1 Learning Goals

After this activity, you should be able to do the following:

  • Plot data points on top of a map using the ggmap() function along with ggplot2 functions
  • Create choropleth maps using geom_map()
  • Add points and other ggplot2 features to a map created from geom_map()
  • Understand the basics of creating a map using leaflet, including adding points and choropleths to a base map


3 Public Electric Vehicle Charging Stations in the United States

As the number of electric vehicles used increases, so do the number of electric vehicle charging stations and ports. The Open Charge Point Interface (OCPI) protocol is a standard for classifying stations according to (1) location, (2) number of charging ports, and (3) number of connectors/plugs (multiple connectors can be on a port, but only one can be used at any given time). Charging equipment is separated into Level 1 charging, Level 2 charging, and DC fast charging, in increasing speed of charging. Read more about these standards here.

According to the Alternative Fuels Data Center (AFDC), as of January 2024, there are 61,124 public electric vehicle charging stations in the United States, and they have 163,191 ports. These numbers have approximately doubled in the last three years:

Image source: Alternative Fuels Data Center.

Where are these stations? How is there distribution across the country related to land area? population? electric vehicle ownership? demographic measures such as socioeconomic, race, and ethnicity?

You can read more about how AFDC collects its data here. Let’s load and display the relevent part of the data:


Exercise 3.1 (Warmup: Stations by state)

  1. Make a chart with the states on the horizontal axis and vertical bars corresponding to the number of stations in each state. Sort the bars so that the state with the most stations is on the left. Add a title.

  2. Now flip the coordinates so that the bars are horizontal and sort them so that the state with the most stations is on top. Then use the options in theme() to align your plot title all the way to the left of the graphic, turn off the background, turn off the grid lines, turn off the tick marks for both axes, and turn off the axis labels for the numeric horizontal axis.

  3. Which version of the chart do you like better? This video on The Glamour of Graphics (you can probably start at 3:45) by William Chase has a lot of good advice about layout, color, and typography for your plots.

ggplot(ev_stations, aes(x=fct_infreq(state)))+
  geom_bar()+
  labs(x="State", title="Number of EV Charing Stations")

ggplot(ev_stations, aes(x=fct_rev(fct_infreq(state))))+
  geom_bar()+
  labs(x="State", title="Number of EV Charing Stations")+
  coord_flip()+
  theme(panel.background=element_blank(),
        axis.ticks=element_blank(),
        axis.title.x=element_blank(),
        axis.text.x=element_blank())+
  geom_text(aes(label=..count..), stat = "count", hjust=0)

Solution.

  1. The second verious is nicer, but took more work. Having the bars go across allows you to see the labels better.


4 Plotting the Charging Station Locations on a Map

We can start by using ggplot2 to make a point plot of all of the stations:

ggplot(data=ev_stations) +
  geom_point(aes(x = longitude, y = latitude), 
             alpha = 0.2, 
             size = .2)

4.1 Example: Superimposing Point Data on Top of a Map

To highlight the geographical nature of this scatterplot, we can superimpose the points on top of a map, using the ggmap() function from the ggmap library. We’ll use the stamen maps from Stadia Maps.

Exercise 4.1 (Sign up for a Stadia API key) To use the Stadia Maps, you’ll each need your own registered API key. To get that, follow the instructions here. Specifically:

  1. Go to client dashboard and sign up for a new account. Verify it by clicking on the email they send you.

  2. Click “manage properties” and then “add an API key.” You can name it something like “RStudio”.

  3. Copy the API key into the code chunk below and delete the eval=FALSE:

register_stadiamaps("5081c024-5bc9-46e4-b756-897ac8e430cd", write = TRUE)

Note: If you prefer Google maps, they also work with ggmap(). We used to be able to easily bring in Google maps, but now you need to have a registered API key for that as well, and they require you to give a credit card number. If you want to do that, see google_key in the help. Then, see the documentation for get_map(). An account does come with a significant amount of free use, but for the required class materials, we are not going to do anything that requires a credit card.

Now we are ready to make an example and walk through the syntax.

# Get the map information
mymap <- get_stadiamap(
    bbox = c(left = -170, bottom = 10, right = -60, top = 66), 
    maptype = "stamen_toner_lite",
    zoom =4)
# Plot the points on the map
ggmap(mymap) + # creates the map "background"
  geom_point(data = ev_stations, 
             aes(x = longitude, y = latitude), 
             alpha = .05, 
             size = .15,
             color="blue")+
  theme_map()

Next, we will walk through the get_get_stadiamap() function arguments. The code below is what was used to get the world map information.

get_stadiamap(
    bbox = c(left = -170, bottom = 10, right = -60, top = 66), 
    maptype = "stamen_toner_lite",
    zoom =4)

bbox

The bbox argument tells it the minimum and maximum latitude and longitude points. So, left is the minimum longitude, right is the maximum longitude, bottom is the minimum latitude, and top is the maximum latitude. One helpful trick is to go to openstreetmap: zoom in on the area of interest, click export, and you will see all the values you need. You may have to modify them slightly, which you can do after your initial plot.

maptype

The maptype tells it the style of the map. Check out the different options by looking in the get_stadiamap help (type ?get_stadiamap in the console).

zoom

When you make a large area, you need to decrease the zoom, otherwise it will take too long to load. So, it’s a good idea to start with a small zoom and you can always make it bigger if you want. This might seem counter-intuitive at first. Think of the zoom level as the level of detail. So, smaller numbers show less detail and larger numbers more detail. A good trick is to go to the stamenmaps webpage and search for the location you are mapping. Then, in the URL, you can see the zoom number. For example, this link is a map of the Olin area: https://maps.stamen.com/#toner/12/37.7706/-122.3782. Notice the number 12 next to /#toner/. That means it is zoomed in at 12.

ggmap()

We save the the map information from get_stadiamap() to a named value and then use it in ggmap():

ggmap(mymap)

The ggmap() function will print the “background” map. Think of it as the providing the canvas on which we will plot. This takes the place of our usual ggplot().

After that, we can use the geom_XXX() functions from ggplot2 that we are used to in order to put points, lines, etc. on top of the map. But, we need to remember to also provide the data we are using in the geom_XXX() function(s) we use, since we do not have the ggplot() function in which to provide it.

# Plot the points on the map
ggmap(mymap) + # creates the map "background"
  geom_point(data = ev_stations, 
             aes(x = longitude, y = latitude), 
             alpha = .05, 
             size = .15,
             color="blue")+
  theme_map()

theme_map()

The last line of the code is theme_map(). This is an optional line from the ggthemes package, but it often makes it look nice.

qmplot()

An alternative is the function qmplot(), which allows you to not specify the bounding box:

qmplot(longitude, latitude, data = ev_stations, maptype = "alidade_smooth", zoom=5,color = I("blue"))


4.2 Additional Resources


4.3 Practice

Exercise 4.2

  1. Add an aesthetic to the example map above that sets the size of the points according to the total_ports variable. What, if anything, can you deduce from this visualization?
  2. Construct a new map of electric vehicle locations in the area around Olin (within a few towns).
  3. In the Olin plot, play with the zoom number. What does it do? (just describe what it does - don’t actually include more than one map).
  4. Try a couple different map types (see get_stadiamap() in help and look at maptype). Include a map with one of the other map types.
  5. Add a point to the map that indicates Olin College and label it appropriately. There are many ways you can do this, but it may be easiest with the annotate() function (see ggplot2 cheatsheet).
# Plot the points on the map
ggmap(mymap) + # creates the map "background"
  geom_point(data = ev_stations, 
             aes(x = longitude, y = latitude, size = total_ports), 
             alpha = .05,
             color="blue")+
  scale_size_continuous(range = c(0.1, 3))+
  theme_map()

# Plot the points on the map
olinmap <- get_stadiamap(
    bbox = c(left = -71.4, bottom = 42.2, right = -71.1, top = 42.35), 
    maptype = "stamen_toner_lite",
    zoom =12)
ggmap(olinmap) + # creates the map "background"
  geom_point(data = ev_stations, 
             aes(x = longitude, y = latitude, size = total_ports), 
             alpha = .5,
             color="blue")+
  scale_size_continuous(range = c(1, 6))+
  theme_map()

# Plot the points on the map
olinmap2 <- get_stadiamap(
    bbox = c(left = -71.4, bottom = 42.2, right = -71.1, top = 42.35), 
    maptype = "stamen_terrain",
    zoom =12)
ggmap(olinmap2) + # creates the map "background"
  geom_point(data = ev_stations, 
             aes(x = longitude, y = latitude, size = total_ports), 
             alpha = .5,
             color="blue")+
  scale_size_continuous(range = c(1, 6))+
  annotate("point", x=-71.2525144, y=42.299611, shape=18, size=3, color="orange")+
  annotate("text", x=-71.2525144, y=42.305, size=3, label="Olin College of Engineering", color="black")+
  theme_map()

Solution.

  1. Can’t tell too much because the map is very busy, but it seems like dense areas tend to have more ports per station, which makes sense
  2. Increases level of detail in the given area


5 Contour Maps

The geom_density_2d and stat_density_2d functions are great for plotting distributions over spatial regions. Here is an example that shows the densities of public electric vehicle charging stations in the United States.

us48 <- get_stadiamap(
    bbox = c(left = -132, bottom = 20, right = -65, top = 55), 
    maptype = "stamen_terrain",
    zoom = 4)

ggmap(us48)+
  geom_density_2d(data=ev_stations, aes(x=longitude,y=latitude),size=.3)+
  stat_density_2d(data = ev_stations, 
    aes(x=longitude,y=latitude, fill = ..level.., alpha = ..level..), 
    size = 0.1, bins = 20, geom = "polygon") + 
  scale_alpha(guide = FALSE) +
  scale_fill_gradient(low = "yellow", high = "red")+
  theme_map()


6 Choropleths

Geographical data needn’t be expressed by latitude and longitude. For choropleth maps, instead of visualizing our data as points with different aesthetics (size, color, transparency, etc.), we color different regions of the maps based on data values. To do this we need to specify both the geometric regions on which the data resides (counties, states, zip codes, etc.), and then wrangle the data so that there is one value per region.

Let’s return to the electric vehicle charging stations and make some choropleths.

6.1 Method 1: geom_map()

The main method we’ll examine is to use geom_map() to create a choropleth map.

6.1.1 Example 1: State Data

The states_map dataset gives information about creating the borders of the US states. The data is retrieved using the map_data() function. Run ?map_data in the console to see more information about what other maps are available. There are also other packages that provide different types of maps.

Let’s take a look at the map and we’ll go through the details afterwards.

# load the data aggregated by state
ev_stations_by_state<-read_csv("http://faculty.olin.edu/dshuman/DS/ev_stations_by_state.csv")
#US states map information - coordinates used to draw borders
states_map <- map_data("state")

# map that colors state by number of charging stations
ggplot(ev_stations_by_state) +
  geom_map(map = states_map,
           aes(map_id = state_name,
               fill = num_stations)) +
  #This assures the map looks decently nice:
  expand_limits(x = states_map$long, y = states_map$lat) + 
  theme_map()

or with a different color scheme:

ggplot(ev_stations_by_state) +
  geom_map(map = states_map,
           aes(map_id = state_name,
               fill = num_stations)) +
  expand_limits(x = states_map$long, y = states_map$lat) + 
  scale_fill_gradientn(colours = c("blue", "white", "red"))+
  theme_map()

Now, let’s look more closely at what each piece of the code below is doing.

ggplot(ev_stations_by_state) +
  geom_map(map = states_map,
           aes(map_id = state_name,
               fill = num_stations)) +
  expand_limits(x = states_map$long, y = states_map$lat) + 
  theme_map()

map: choose a map

The map argument tells R at which level to create the map. Really, it tells it how to draw all the borders This is a very special data set. According to the geom_map() documentation, it is a “data frame that contains the map coordinates … It must contain columns x or long, y or lat, and region or id.” We are using the map_data() function to create the map file (see above for more detail). You can open the map data, states_map, and see that it adheres to the rules.

map_id: connect map ID/region variable to data being plotted

The map_id inside of aes() is a required aesthetic for the geom_map() geom. It tells R which variable is the region/id variable, in this case the state. It connects the region or id from the map (region variable in states_map dataset, in this example) to the dataset being plotted (state_name in ev_stations_by_state, in this example). So state_name needs to have the same form as region.

use ggplot2 features

We tell it to fill in the states by the variable num_stations, the number of electric vehicle charging stations in each state. With the geom_map() geom, it will fill in the borders of the regions we defined in the map argument.

expand_limits()

Use expand_limits() to assure that the map covers the entire area it’s supposed to. We put the longitude variable from states_map for the x argument and the latitude variable from states_map for the y argument to assure the map stretches across the entire range of longitudes and latitudes in the map.

theme_map()

This is a personal preference, but theme_map() often makes the map look nicer.

add ggplot2 layers

You can add any of the ggplot2 layers on top of this map. For example, we can add a title and the station locations as points:

ev_stations48<-filter(ev_stations,!state %in% c("PR","AK","HI") )
ggplot(ev_stations_by_state) +
  geom_map(map = states_map,
           aes(map_id = state_name,
               fill = num_stations)) +
    geom_point(data = ev_stations48,
           aes(x = longitude, y = latitude),
           size = .05,
           alpha = .2, 
           color = "goldenrod")+
  lims(x=c(-125,-65),y=c(23,50))+
  theme_map()+
  ggtitle("Public Electric Vehicle Charging Stations in the United States")

6.1.2 Example 2: County Data

This data set did not come with the county name or FIPS code for each station. So, using the process outlined here, I had R associate a county with the latitude-longitude coordinates for each stations, tacked the county name and FIPS code onto the data set, and tallied up the stations by county:

ev_stations_by_county<-read_csv("http://faculty.olin.edu/dshuman/DS/ev_stations_by_county.csv")

Now we can load a county map and make a choropleth:

county_map<-socviz::county_map # from socviz library
ggplot(ev_stations_by_county)+
  geom_map(data=ev_stations_by_county,map=county_map,aes(map_id=fips, fill=num_stations))+
  expand_limits(x = county_map$long, y = county_map$lat)+
  theme_map()+
  theme(legend.position = "right")+
  labs(fill=element_blank(),title="Number of Public Electric Vehicle Charging Stations by United States County")

With this coloring scale, it is nearly impossible to tell a county that has 0 stations from a county that has 10 stations or a county that has 100 stations. One way to address is this is to make a new discrete (categorical) variable by using cut() to group the counties into bins according to their number of stations:

# define appropriate (& nicely labeled) breaks
ev_stations_by_county$num_stations_cat <- cut(ev_stations_by_county$num_stations, 
                      breaks=c(-1,.99,10,100,seq(500,4000,by=500)), 
                      labels=c("0", "1-10", "11-100", "101-500", 
                               "501-1000", "1001-1500", "1501-2000", "2001-2500", "2501-3000", "3001-3500","3501-4000"),
                      include.lowest=TRUE)

Now we can make a second map that is easier to read:

ggplot(ev_stations_by_county)+
  geom_map(data=ev_stations_by_county,map=county_map,aes(map_id=fips, fill=num_stations_cat))+
 scale_fill_manual(values = rev(brewer.pal(8,"RdBu")))+
  expand_limits(x = county_map$long, y = county_map$lat)+
  theme_map()+
  theme(legend.position = "right")+
  labs(fill=element_blank(),title="Number of Public Electric Vehicle Charging Stations by United States County")

6.1.3 Additional Resources


6.2 Alternative Choropleth Methods

There are plenty of other methods available to make choropleths in R.

Let’s demonstrate just three additional methods.

6.2.1 Alternative 1: ggplot + geom_polygon

mapping_data <- ev_stations_by_county %>% 
  mutate(id = fips) %>% 
  left_join(county_map, by = "id")
ggplot(mapping_data, aes(x = long, y = lat, fill = num_stations_cat, group = group)) + 
  coord_equal() + 
  geom_polygon(color = NA)+
  theme_map()+
  theme(legend.position = "right")+
  scale_fill_manual(values = rev(brewer.pal(8,"RdBu")))+
  labs(fill=element_blank(),title="Number of Public Electric Vehicle Charging Stations by United States County")

6.2.2 Alternative 2: plot_usmap

# This function is in the usmap package
plot_usmap(data = ev_stations_by_county, values = "num_stations_cat",color=NA)+
  scale_fill_manual(values = rev(brewer.pal(8,"RdBu")))+
  labs(fill=element_blank(),title="Number of Public Electric Vehicle Charging Stations by United States County")+
  theme(legend.position = "right")

6.2.3 Alternative 3: choroplethr Package

The county_choropleth function in the choroplethr package requires the variable of interest to be stored as value and the FIPS code to be stored in region. The following code does this for the charging station data:

ev_stations_by_county <- mutate(ev_stations_by_county, value=num_stations_cat,region=as.double(fips))
county_choropleth(ev_stations_by_county)+
scale_fill_manual(values = rev(brewer.pal(8,"RdBu")))+ labs(fill=element_blank(),title="Number of Public Electric Vehicle Charging Stations by United States County")


We can also zoom in on specific states:

county_choropleth(ev_stations_by_county, state_zoom="massachusetts")+
 scale_fill_manual(values = rev(brewer.pal(8,"RdBu")))+ labs(fill=element_blank(),title="Number of Public Electric Vehicle Charging Stations by Massachusetts County")

For this one , we can try increasing the color resolution by giving it the continuous variable directly and letting the package choose the discrete colors:

ev_stations_by_county2 <- mutate(ev_stations_by_county, value=num_stations)
county_choropleth(ev_stations_by_county2, state_zoom="massachusetts")+ labs(fill=element_blank(),title="Number of Public Electric Vehicle Charging Stations by Massachusetts County")


6.3 Practice


The examples above did not account for population of each state in the map. The ev_stations_with_2018_pop_est dataset includes a new variable, ev_stations_per_10000, that gives the number of public electric vehicle charging stations per 10,000 people.

Exercise 6.1 Create a choropleth map that shows the number of electric vehicle charging stations per 10,000 people on a map of the US. Use a new fill color, add an informative title for the plot, and include a caption that says who created the plot (you!). Make a conclusion about what you observe.

ggplot(ev_stations_with_2018_pop_est) +
  geom_map(map = states_map,
           aes(map_id = state_name,
               fill = ev_stations_per_10000)) +
  #This assures the map looks decently nice:
  expand_limits(x = states_map$long, y = states_map$lat) +
  labs(title="EV Charging Stations Per 10k People", caption = "Created by Ian Lum")+
  theme_map()

::: {.solution} We see that California had such a larger number of stations mainly because it is a very populous state. However, Vermont actually appears to have the most EV stations after adjusting for population. :::


7 Dynamnic Maps with leaflet

Leaflet is an open-source JavaScript library for creating maps. It can be used outside of R, but we will only discuss using the leaflet library in R.

This library uses a different plotting framework from ggplot2 although it still has a tidyverse feel due to its use of the pipe, %>% and the way it adds layers to the plot, just like in ggplot2.

7.1 Steps to Create a Map

  1. Create a map widget by calling leaflet() and telling it the data to use.
  2. Add a base map using addTiles() (the default) or addProviderTiles().
  3. Add layers to the map by using layer functions (e.g. , addMarkers(), addPolygons()) to modify the map widget.
  4. Repeat step 3 as desired.
  5. Print the map widget to display it.

7.2 Creating a Map with Markers/Points

Below, we create a basic map and add points of interest. You can use a “two-finger scroll” to zoom in and out. WARNING: DO NOT FORGET THE ~ BEFORE THE VARIABLE NAMES!!!

ev_stations_needham_wellesley<-filter(ev_stations,state=="MA",city %in% c("Needham","Wellesley")) # you can add more places, but it may take a while to run
leaflet(data = ev_stations_needham_wellesley) %>% #base plot
  addTiles() %>% #base map - default is openstreet map 
  addMarkers(lng = ~longitude, 
             lat = ~latitude, 
             label = ~station_name) 

We can change just about everything about our map:

  • The base map with addProviderTiles().

    • To see all available provider base maps, type providers in the console.

    • To access those maps, use providers$PROVIDERNAME inside the addProviderTiles() function, where PROVIDERNAME is one of those listed providers. When you type provider$ a list should show up that you can click on.

  • The marker type by using addCircles() instead of addMarkers() - Search addControl in the Help or type ?addControl into the console to see what all the arguments mean and how you can change them.

Here is one more example showing some of these options:

leaflet(data = ev_stations_needham_wellesley) %>%
  addProviderTiles(providers$CartoDB.Positron)%>%
  addCircles(lng = ~longitude,
             lat = ~latitude, 
             label = ~station_name, 
             weight = 10, 
             opacity = 1, 
             color = col2hex("darkblue")) 

Some other notes:

  • All variable arguments are preceded with a tilde, ~.

  • The weight argument tells it how thick to make the lines or points, pixels.

  • The opacity argument is the transparency, like the alpha argument in ggplot2.

  • Colors need to be in “hex” form. The col2hex() function from the gplot library will do that. The colors also need to be valid R colors.

7.3 Choropleth Layers with addPolygons()

For making maps with borders (like choropleth maps), the functions can receive the base maps with spatial data a few different ways. In the example here, we use functions from the sf package to get the data in the right form. In the demo video listed under Resources below, a data.frame is translated to an sf object. Hopefully those two sets of instructions fit most of the cases you are interested in.

In the data code chunk below, the function st_read() downloads the shape file for the counties of North Carolina, which is included in the sf package.

# North Carolina births and sids deaths 
nc <- st_read(system.file("shape/nc.shp", package="sf"), 
              quiet = TRUE) %>% 
  mutate(sid_per_1000birth_79 = SID79/BIR79*1000)

# Compute centroid (center) of each county
county_centroid_lat_long <- as_tibble(st_coordinates(st_centroid(nc))) %>% 
  rename(centr_long = X, 
         centr_lat = Y)

# County names and sid_per_1000birth_79
nc_centroid_county <- st_drop_geometry(nc) %>% 
  select(NAME, sid_per_1000birth_79)

# All centroid level info
nc_centroid <-  county_centroid_lat_long %>% 
  bind_cols(nc_centroid_county)

The dataset has number of births and number of SIDS cases in each county of North Carolina from 1974-1979 and 1979-1984. We computed a variable called sid_per_1000birth_79 which is the number of SIDS cases per 1000 births in 1979.

Below, the the NAMES and geometry variables from the first five rows of the data are printed out. The geometry variable contains information about how to plot the boundaries of the counties. Open the dataset and examine the geometry variable. The leaflet function knows that the geometry variable contains this special information.

nc %>% 
  select(NAME, geometry) %>% 
  slice(1:5)
## Simple feature collection with 5 features and 1 field
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -81.74107 ymin: 36.07282 xmax: -75.77316 ymax: 36.58965
## Geodetic CRS:  NAD27
##          NAME                       geometry
## 1        Ashe MULTIPOLYGON (((-81.47276 3...
## 2   Alleghany MULTIPOLYGON (((-81.23989 3...
## 3       Surry MULTIPOLYGON (((-80.45634 3...
## 4   Currituck MULTIPOLYGON (((-76.00897 3...
## 5 Northampton MULTIPOLYGON (((-77.21767 3...

To learn more about the data type nc into the Help menu or ?nc into the console.

This plot puts the map with North Carolina county borders on top of the Open Street Map. This map isn’t that interesting by itself.

leaflet(nc) %>% 
  addTiles() %>% 
  addPolygons()

Now, let’s use some of the data to enhance the graph by creating a choropleth map where we color the counties by sid_per_1000birth_79. In ggplot(), we can map a variable to color or fill inside the aesthetic, but in the leaflet functions we cannot do this. Instead, we need to create a variable of hex color names that tells it how to color or fill. Thankfully, there are functions that help us do that! Let’s walk through the detail of how we created the graph below.

#creates a function that can be used to map a variable to
#the color palette using viridis
pal <- colorNumeric("viridis", 
                     domain = nc$sid_per_1000birth_79) 

leaflet(nc) %>% 
  addTiles() %>% 
  addPolygons(
    fillColor = ~pal(sid_per_1000birth_79), #fills according to that variable
    fillOpacity = 0.7) #like alpha in ggplot

The colorNumeric() function returns a function that maps a variable’s values to colors in the given palette, in this case “viridis”. So, pal() is a function. We can then use that function inside addPolygons(). If we apply the function to the sid_per_1000birth_79 variable, it returns a variable of hex colors and the variable of colors is used to fill the counties.

Below we print out what happens when pal() is applied to sid_per_1000birth_79. We can see that it returns hex colors. So, the variable in the fillColor argument inside addPolygons() above, is a variable of hex color names.

head(pal(nc$sid_per_1000birth_79))
## [1] "#440154" "#C0DF25" "#38588C" "#2A778E" "#34618D" "#26828E"

Again, this is different from ggplot(), where we could map a variable to color or fill and it would do the translating of variable to color scale for us. In the leaflet functions, we have to explicitly provide the colors in a variable. The colorNumeric() command helps you do that. colorBin(), colorQuantile(), and colorFactor() are other functions you might need to use depending on the type of variable you are trying to map to colors.

There are many customizeable options in leaflet, much like ggplot(). Here is a commented set of code to point out some useful functions and arguments:

#creates a function that can be used to map a variable to
#the color palette using viridis
pal <- colorNumeric("viridis", 
                     domain = nc$sid_per_1000birth_79) 

leaflet(nc) %>% 
  addTiles() %>% 
  addPolygons(
    #skips drawing the borders:
    stroke = FALSE, 
    #fills according to variable of hex colors:
    fillColor = ~pal(sid_per_1000birth_79), 
    #changes transparency, like alpha in ggplot
    fillOpacity = 0.7, 
    #how much to simplify the plot when zooming:
    smoothFactor = 0.5, 
    #changes what happens to the shape when we mouse over it
    highlight = highlightOptions(weight = 5, 
                                 color = "black",
                                 fillOpacity = 0.9,
                                 bringToFront = FALSE)) %>%
  addCircles(data = nc_centroid,
             lng = ~centr_long, 
             lat =  ~centr_lat, 
             # label that appears when you click on the marker,
             # in this case county name and sid_per_1000birth_79
             # rounded to the 2nd decimal
             popup = ~paste(NAME,": ",
                            round(sid_per_1000birth_79, 2),
                            sep=""),
             radius = 2) %>% 
  
  # Add a legend
  addLegend(pal = pal, 
            values = ~sid_per_1000birth_79, 
            opacity = 0.5, 
            title = NULL,
            position = "bottomright") 

7.4 Additional Resources

7.5 Practice

Exercise 7.1 In this exercise, you are going to create a single map of some of your favorite places! The end result will be one map.

  1. Using the tibble() function, create a data set that has 10-15 rows of your favorite places. The columns will be the name of the location, the latitude, the longitude, and a column that indicates if it is in your top 3 favorite locations or not.

  2. Create a map that uses circles to indicate your favorite places. Label them with the name of the place. Choose the base map you like best. Color your 3 favorite places differently than the ones that are not in your top 3. Add a legend that explains what the colors mean.

  3. If there are other variables you want to add that could enhance your plot, do that now. You could also connect all your locations together with a line in a meaningful way (you may need to order them differently in the original data).

places <- tribble(
  ~name, ~lat, ~long, ~top3,
  "Olin", 42.299611, -71.2525144, FALSE,
  "Bocado", 42.2958182, -71.297112, TRUE,
  "Laughing Monk", 42.295788667666955, -71.29760080747083, FALSE,
  "French Press", 42.2819488875885, -71.23656883286988, TRUE,
  "Trader Joe's", 42.292816108312145, -71.23465745028308, FALSE,
  "Cafe Fresh Bagel", 42.29357541368595, -71.23513877879792, FALSE,
  "Little Big Diner", 42.33097980148705, -71.19370982227535, FALSE,
  "Dave's Hot Chicken", 42.31360499068687, -71.21296140476686, FALSE,
  "Penang", 42.36687603695944, -71.23797019480519, TRUE,
  "The Local", 42.307582797592325, -71.27865511599475, FALSE
)
placemap <- get_stadiamap(
    bbox = c(left = -71.34, bottom = 42.26, right = -71.15, top = 42.4), 
    maptype = "stamen_toner_lite",
    zoom =12)
ggmap(placemap,
      base_layer = ggplot(
        places,
        aes(x = long, y = lat, color = top3, label = name)
      )) + # creates the map "background"
  geom_point(alpha = .7, size=3)+
  geom_text(vjust=-1)+
  labs(title="Ian's Favorite Places", color="Top 3 Places")+
  theme_map()